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ABSTRACT 

The  sensitivity  of  the  onset  and  the  location  of  vortex  breakdowns  in  concentrated  vor¬ 
tex  cores,  and  the  pronounced  tendency  of  the  breakdowns  to  migrate  upstream  have  been 
characteristic  observations  of  experimental  investigations;  they  have  also  been  features  of 
numerical  simulations  and  led  to  questions  about  the  validity  of  these  simulations.  This 
behavior  seems  to  be  inconsistent  with  the  strong  time-like  axial  evolution  of  the  flow,  as  ex¬ 
pressed  explicitly,  for  example,  by  the  quasi-cylindrical  approximate  equations  for  this  flow. 
An  order-of-magnitude  analysis  of  the  equations  of  motion  near  breakdown  leads  to  a  mod¬ 
ified  set  of  governing  equations,  analysis  of  which  demonstrates  that  the  interplay  between 
radial  inertial,  pressure,  and  viscous  forces  gives  an  elliptic  character  to  these  concentrated 
swirling  flows.  Analytical,  asymptotic,  and  numerical  solutions  of  a  simplified  non-linear 
equation  are  presented;  these  qualitatively  exhibit  the  features  of  vortex  onset  and  location 
noted  above. 
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1.  Introduction 


The  study  of  the  breakdown  of  the  cores  of  strong  longitudinal  vortices  has  been  of 
recurring  interest  for  many  years.  Breakdown  may  be  described  as  the  abrupt  change  in  the 
structure  of  the  core  of  a  vortex,  often  characterized  by  the  presence  of  a  free  stagnation 
point  on  the  axis  of  the  vortex  and  a  corresponding  divergence  of  the  stream  surfaces  near  the 
axis,  and  in  some  cases  a  region  of  reversed  flow.  In  addition  to  being  of  theoretical  interest, 
vortex  breakdown  has  important  technological  applications,  both  aerodynamic  and  non- 
aerodynamic  [1].  In  some  applications  vortex  breakdown  is  beneficial,  for  example,  in  leading 
to  enhanced  mixing  or  momentum  or  heat  exchange;  or  destructive,  as  in  degradation  of 
aerodynamic  performance.  Vortex  breakdown  may  play  an  important  role  in  vortex  dynamics 
in  general,  and  perhaps  also  in  transition  to  turbulence.  Given  the  importance  of  vortex 
breakdown,  it  is  quite  remarkable  that  the  basic  underlying  mechanism  leading  to  breakdown 
is  not  yet  known  or  generally  accepted!  This  is  not  for  lack  of  postulated  explanations  the 
literature  abounds  with  them.  More  remarkable  yet  is  that  these  explanations  are  so  vastly 
different  in  what  they  postulate  to  be  the  underlying  mechanism;  these,  variously,  affirm 
(their  principal  claimants  indicated  parenthetically)  that: 

(i)  vortex  breakdown  is  similar  to  the  separation  of  a  two-dimensional  boundary  layer 
(Hall  [2,3],  Mager  [4]); 

(ii)  vortex  breakdown  is  a  consequence  of  hydrodynamic  instability  (Ludwieg  [5,6]); 

(iii)  vortex  breakdown  is  dependent  on  the  existence  of  a  critical  state,  and  is  a  finite 
transition  between  states,  analogous  to  a  hydraulic  jump  (Squire  [7],  Benjamin  [8,9]); 

(iv)  vortex  breakdown  is  like  a  solitary  wave,  or  soliton,  the  result  of  the  trapping  of  long, 
weakly  nonlinear  waves  propagating  in  nearly  critical  swirling  flows  (Leibovich  [10,11]). 

To  a  greater  or  lesser  degree  each  of  these  theories  fails  to  explain  fully,  or  adequately, 
breakdown  in  all  its  aspects.  Experiments  and  numerical  simulations  have  failed  to  resolve 
the  matter  satisfactorily.  How  can  one  account  for  this  failure?  A  number  of  reasons  present 
themselves.  First,  experiments  and  measurements  are  difficult.  There  is  a  tendency  for  the 
breakdown  to  migrate  back  and  forth  in  the  test  section  or  on  the  aerodynamic  surface.  Mea¬ 
surements,  whether  invasive  (e.g.  hot  wires)  or  non-invasive  (e.g.  LDV),  and  flow  visualiza¬ 
tion  are  difficult  to  obtain  and  to  interpret  in  these  potentially  unsteady,  three-dimensional 
flows.  Over  the  years  some  of  the  seemingly  different  postulated  physical  mechanisms  of 
vortex  breakdown  have  been  shown  to  be  equivalent  or  related,  leading  to  the  same  or  sim¬ 
ilar  criteria  for  breakdown.  To  the  extent  that  the  various  theories  predict  a  criterion  for 
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vortex  breakdown  the  predictions  are  quite  similar  [12],  namely  that  tansup  —  1  v/w  as  57 
(de  or  v/w  ~  1.5,  where  v  and  w  are  the  swirl  and  axial  velocities,  respectively,  so  this 
cannot  be  used  to  distinguish  among  or  between  the  various  physical  mechanisms  and  their 
corresponding  theoretical  formulations. 

As  for  numerical  simulations  of  breakdown  using  the  full  Navier-Stokes,  or  inviscid  Euler 
equations,  until  recently  these  generally  required  the  assumption  of  axisymmetry  and  steadi¬ 
ness  of  the  flow  [13-16].  This  is  not  true  of  the  more  recent  calculations,  for  example,  those 
by  Kuruvila  &  Salas  [17],  Spall  et  al.  [18,19],  Breuer  &  Hanel  [20],  and  Menne  &  Liu  [21],  all 
of  which  solve  the  unsteady,  three-dimensional,  fully  viscous  Navier-Stokes  equations  (but 
not  all  of  them  are  necessarily  time-accurate  solutions).  Common  characteristics  of  many  of 
the  earlier  and  these  most  recent  numerical  simulations,  and  for  that  matter  of  the  experi¬ 
ments  as  well,  are:  (i)  extreme  sensitivity  to  flow  parameters  (such  as  swirl  velocity  ratio, 
external  axial  velocity  variation,  or  pressure  gradient);  (ii)  suddenness  of  breakdown  (i.e.  no 
evidence  of  the  growth  of  an  instability  of  the  basic  swirling  flow);  and  (iii)  a  tendency  for 
the  breakdowns  to  migrate  upstream  to  the  initial  station  or  computational  cell.  In  some 
simulations  [18]  this  is  avoided  by  specifying  the  external  axial  velocity  gradient. 

This  last  characteristic  of  these  computations  is  the  case,  at  least,  for  unbounded  or 
“open”  flows,  e.g.  aerodynamic  flows  and  flows  in  diverging  channels.  Numerically  simulated 
breakdowns,  usually  in  the  form  of  steady  axisymmetric  bubbles,  that  occur  in  bounded  or 
confined  swirling  flows,  for  example,  in  finite  closed  cylinders  with  one  endwall  rotating  (Lugt 
[22,23],  Neitzel  [24])  or  flow  in  the  gap  between  rotating  spheres  (Bar-Yoseph  et  al.  [25])  do 
not  exhibit  this  behavior,  nor  do  experiments  on  these  flows  (Escudier  [26]).  This  tendency 
of  the  breakdowns  to  migrate  to  the  initial  station  has  caused  questions  to  be  raised  about 
the  validity  of  such  numerical  simulations  [27].  Another  potential  concern  or  question  in  the 
simulations  of  breakdown  is  the  appropriate  downstream  boundary  conditions  to  impose  in 
order  to  “close”  the  flow  domain,  a  question  whether  one  uses  the  elliptic  Euler  or  Navier- 
Stokes  equations.  What  makes  these  issues  more  troubling  is  the  question  of  the  extent  to 
which  difficulties  with  the  numerical  simulations,  including  the  slow  rates  of  convergence  to 
final  steady  states  often  encountered,  reflect  physical  reality,  or,  in  other  words,  are  intrinsic 
aspects  of  the  onset  and  location  of  breakdown. 

This  paper  addresses  and  attempts  to  answer  one  of  the  above  questions:  the  extreme 
sensitivity  of  the  breakdown  location  to  certain  flow  parameters  and  the  pronounced  ten¬ 
dency  for  the  breakdown  to  migrate  upstream.  The  theoretical  explanation  for  these  effects  is 
given  in  Section  2  where  the  “quasi-cylindrical  equations”  (Hall  [2])  are  analyzed  in  some  de¬ 
tail.  These  equations  are  simplified  forms  of  the  Navier-Stokes  equations  derived  by  making 
boundary  layer-like  assumptions  and  approximations,  and  like  those  equations  are  parabolic, 
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reflecting  a  strong  “time-like”  character  to  the  spatial  evolution  of  these  thin  strongly-swirling 
vortex  cores.  We  show  that  near  breakdown,  in  the  radial  momentum  equation,  which  is  a 
purely  inviscid  equation  in  the  quasi-cylindrical  approximation,  the  radial  viscous  diffusion 
term  must  be  retained.  An  analysis  of  this  modified  set  of  equations,  which  would  seem,  like 
the  original  equations,  to  be  parabolic  since  no  axial  viscous  diffusion  terms  appear,  shows 
that  in  fact  they  have  elliptic-like  characteristics.  This  is  explicitly  shown  by  demonstrating 
that  the  equations  may  be  recast  in  the  form  of  an  elliptic-integro-partial  differential  equa¬ 
tion  for  the  axial  velocity  variation.  In  Section  3  a  simplified,  but  still  non-linear,  equation 
is  employed  as  a  model  of  the  full  integro-differential  equation  and  its  full  solution  given. 
Included  in  this  section  is  also  a  matched  asymptotic  expansion  of  the  equation  in  the  limit 
of  large  Reynolds  numbers.  (This  solution  yields  simpler  expressions  for  the  axial  velocity 
variation  near  breakdown  and  its  explicit  dependence  on  flow  and  fluid  parameters.)  This 
section  concludes  with  an  order-of-magnitude  analysis  which  explicitly  relates  the  constant 
appearing  in  the  simplified  equation  to  the  primary  flow  parameters:  the  core  Reynolds 
number,  the  flow  divergence,  the  swirl  velocity  ratio,  and  the  (non-dimensional)  core  circu¬ 
lation.  In  Section  4  explicit  results  of  the  exact  analytic  solution  presented  in  Section  3  for 
a  representative  range  of  values  of  the  flow  parameters  are  presented  and  discussed. 

2.  Mathematical  Model 

We  assume  the  flow  is  steady,  axisymmetric  and  laminar.  The  governing  Navier-Stokes 
equations  in  a  cylindrical  coordinate  system  r,  0,  z,  with  corresponding  velocity  components 
u,  u,  w  (the  overbars  denote  dimensional  variables)  are 
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We  define  the  following  non-dimensional  variables 


“=(?)il’  p  = 


P- Poo 

pW$  ’ 


r=S ’ 


Z=L' 


where  L  is  a  characteristic  scale  in  the  z-direction,  8  is  a  characteristic  core  radius,  Wo  a 
characteristic  axial  flow  velocity  (e.g.,  free-stream  velocity  at  the  edge  of  the  core),  the 
uniform  static  pressure  far  from  the  vortex,  p  the  density,  and  v  the  kinematic  viscosity.  In 
terms  of  these  quantities  the  governing  equations  become 


13.  .  dw 
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dr 
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dv  uv  1  fL\\d(ld.  A  / 8\2  d2vl 
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dw  dw  dp  1  (L\  1  d  .  dw  ( 8\2  d2w 

u!b+w!h  ~  +  \l)  +  (z)  ~d^y  ^ 

where  Res  =  WqS/u  is  the  Reynolds  number  based  on  the  characteristic  core  radius. 

Equations  (2.6)  -  (2.9)  contain  two  parameters:  8fL  and  Res.  At  large  distances  upstream 
of  breakdown  the  vortex  develops  on  a  length  scale  L  which  is  much  larger  than  8,  so  8/ L  <C  1 
in  this  region.  The  core  Reynolds  numbers  for  which  vortex  breakdowns  occur  are  large, 
but  only  modestly  so,  of  the  order  of  hundreds  ([14-16],  [18]).  We  see  immediately  that 
if  8/L  =  c  <C  1  that  the  axial  viscous  diffusion  terms  in  the  radial,  azimuthal,  and  axial 
momentum  equations,  (2.7)  -  (2.9),  are  negligible  compared  to  the  other  terms  in  these 
equations.  If  we  take  the  formal  limit  c  — ►  0,  with  Res  held  fixed,  and  without  assuming 
that  eRes  is  either  large  or  small,  then  (2.6)  -  (2.9)  reduce  to 


13.  .  dw 

;a;(ru)  +  Tz =  °- 


dv  dv  uv  1  f  L\  \  d  (\  d  ,  A] 

Ulh  +  W~dl  +  T~  Ws  (?)  *  (;*(n’)jJ  ’ 


(2.10) 


(2.11) 

(2.12) 
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(2.13) 


dw  dw  dp  1  ( L\  [l  d  (  dw\ 

U  dr  +  W  dz  dz+  Res  \  8 )  [r  dr  \  dr  J 

Equations  (2.10)  -  (2.13)  are  the  so-called  quasi-cylindrical  equations,  derived  and  extensively 
analyzed  by  Hall  [3].  In  contrast  to  Eqs.  (2.6)  -  (2.9)  which  are  elliptic,  these  quasi-cylindrical 
equations  are  parabolic  equations,  which  can  be  solved  by  time-marching  schemes.  They 
differ  from  the  ordinary  boundary  layer  equations  only  by  the  inviscid  radial  momentum 
equation,  (2.11),  which  contributes  a  non-zero  normal  pressure  gradient.  These  equations 
were  solved  numerically  by  Hall  [27]  to  predict  the  location  of  the  onset  of  breakdown. 
Because  upstream  influence  has  been  neglected  they  cannot  be  expected  to  predict  the  details 
of  the  breakdown  itself. 

As  we  approach  breakdown  the  scale  of  2-derivatives  becomes  comparable  to  that  of  r- 
derivatives,  so  8/L ,  even  if  small,  is  not  expected  to  be  infinitesimal.  Hence,  returning  to 
Eqs.  (2.6)  -  (2.9),  for  this  near-breakdown  region  we  again  drop  all  terms  of  0[(6/L)2],  but 
keep  in  addition  to  the  ( L/8)ReJ 1  terms  also  terms  of  0(8/L).  (In  other  words,  we  are 
assuming  that  8/L  is  small,  but  not  infinitesimal,  so  terms  of  0[(8/L )2]  may  be  neglected, 
but  not  those  of  0(8/ L).)  There  is  only  one  such  term  in  Eqs.  (2.6)  -  (2.9),  the  second  term 
on  the  right  hand  side  of  (2.7);  this  equation  then  becomes 


dp  1  (8\  6  (\  d  . 


(2.14) 


The  other  equations  remain  unchanged,  viz.  (2.10),  (2.12)  and  (2.13).  (Eq.  (2.14)  as  a 
replacement  for  (2.11)  in  the  near-breakdown  region  may  be  more  formally  justified  by  an 
appropriate  rescaling  of  variables  or  by  a  multiple-scale  argument.)  In  dimensional  form  the 
new  set  of  reduced  equations  for  the  near-breakdown  region  is  then 


Id,,  dw 
-rTr(ru)  +  &  =  °’ 

(2.15) 

du  du  v2  1  dp  d  ( 1  d  .  .\ 

u»+”’aI-7  =  -?a;  +  "a;(;s;(r“)j’ 

(2.16) 

dv  dv  uv  d  ( 1  1  , 

“3?  +  “fo  +  T  =  "dr  Vr  a?(r”)J  ’ 

(2.17) 

dw  dw  1  dp  v  d  (  dw\ 

U  dr  +  W  dz  p  dz  +  r  dr  \  dr  ) 

(2.18) 

These  equations  differ  from  the  full  axisymmetric  Navier-Stokes  equations  only  by  the  omis¬ 
sion  of  all  the  axial  viscous  diffusion  terms  and  the  inertia  terms  in  the  radial  momentum 
equation.  They  differ  from  the  quasi-cylindrical  equations,  (2.10)  -  (2.13),  only  by  the  inclu¬ 
sion  of  the  radial  viscous  diffusion  term  in  the  radial  momentum  equation.  Since,  as  for  this 


5 


modified  set  of  reduced  equations,  no  axial  viscous  diffusion  terms  appear  in  (2.15)  -  (2.18) 
we  would  expect  this  set  to  be  parabolic  as  well.  As  we  shall  see  below,  however,  this  is  not 


Combining  (2.15)  and  (2.16)  we  can  write 


v2  _  1  dp  d2w 

r  pdr  drdz 


(2.19) 


If  we  integrate  this  with  respect  to  r  from  0  to  r,  and  then  differentiate  the  resulting  expres¬ 
sion  with  respect  to  z ,  it  becomes 

( dp\  ( dp\  r  1  dv2  [d2to  d2w 

Wr=0"  r !hdr  +  ti  [a?|r_  d?lr=0]  •  (2,20) 

Multiplication  of  (2.17)  by  2v,  and  combining  terms,  leads  to  the  expression 


0,11.  d(r2v2)  n  ,  alia.  .  i 

u— (r  v  )  +  w— —  =  2i /r  l-^Ml  , 


8_  l_d_ 

dr  r  dr 


(2.21) 


which  can  be  used  to  eliminate  dv2/dz  in  the  integral  term  in  (2.20),  resulting  in,  dropping 


dp.  dp  r  1  fu\  d  2  2  r  1  /u\  d  [l  d  .  .1  ,  \d2w  d2w 

=  PL  U)  7r{r  v  )dr~2l,h  r  U  S  ba;<r'’)]  * +"  [a?  _  a?H 

(2.22) 

Equation  (2.22)  exhibits  the  various  mechanisms  that  can  lead  to  the  development  of  an 
adverse  pressure  gradient  along  the  core  axis  and  consequent  axial  retardation  of  the  flow. 
The  first  term  on  the  right  hand  side  is  an  inviscid  term  and  was  first  identified  in  this  context 
by  Hall  [3]  from  an  inviscid  analysis  of  the  quasi-cylindrical  equations.  It  is  in  the  nature  of  a 
flow  divergence-swirl  interaction  term,  demonstrating  how,  for  u  >  0,  arising  from  whatever 
cause,  and  positive  dk/dr ,  where  k  =  rv,  which  is  generally  the  case,  the  pressure  gradient 
along  the  axis  can  exceed  that  for  larger  r,  and  in  particular  that  at  the  core  edge.  If  T  and 
R  are  characteristic  measures  of  the  core  circulation  (i.e.,  rv  evaluated  at  the  core  edge)  and 
radius,  and  we  set  u/w  =  a,  a  measure  of  flow  divergence,  an  order-of-magnitude  analysis 
of  this  integral,  evaluated  at  the  core  edge,  yields  (Hall  [28]) 

?E\  _  ?E\  ^  (oov 

dz\r=0  dz\r=R  R3  .  (2.23) 

Thus,  if  a  >  0,  i.e.  the  stream  surfaces  diverge  in  the  core,  then  the  pressure  gradient  on  the 
axis  exceeds  that  along  the  core  edge  by  an  amount  that  can  be  very  large  for  a  concentrated, 
intense  vortex  core. 

The  two  new  terms  on  the  right  hand  side  of  (2.22)  are  both  viscous  terms.  The  first  of 
these,  the  second  one  on  the  right  hand  side,  represents  the  interaction  of  swirl  and  radial 


(2.23) 
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diffusion  of  swirl  momentum,  and  also  acts  to  raise  the  pressure  on  the  axis.  The  last 
term  on  the  right  has  the  form  of  an  axial  diffusion  term,  which  has  arisen  not  from  the 
inclusion  of  axial  diffusion  terms  in  the  governing  equations,  because  there  are  conspicuously 
absent  from  Eqs.  (2.15)  -  (2.18),  but  from  the  continuity  equation.  Thus  while  the  original 
quasi-cylindrical  equations  are  parabolic,  their  modification,  by  the  inclusion  of  the  radial 
diffusion  term  in  the  radial  momentum  equation,  has  the  effect  of  introducing  a  term  which 
has  the  appearance  of  an  axial  diffusion  term!  This  term  gives  an  elliptic  character  to 
the  equations.  More  precisely,  (2.22)  is  an  elliptic-integro-partial  differential  equation.  In 
addition  to  showing  how  upstream  influence  enters  the  vortex  breakdown  problem,  Eq.  (2.22) 
also  exhibits  explicitly  how  flow  divergence,  a  =  u/w,  and  the  swirl  ratio,  v/w,  independently 
affect  the  axial  pressure  gradient. 

We  can  see  more  immediately  the  effect  of  the  terms  on  the  right  hand  side  of  (2.22)  on 
the  axial  velocity  by  using  (2.18)  to  eliminate  dp/dz  from  (2.22),  obtaining 


(2.24) 


0  r  \  ( v\  d  \\  d ,  .  , 

+2i/  /  -  (  — )  - ^-(rv)  dr  -  v 
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Equation  (2.24)  contains  as  unknowns  the  three  velocity  components.  To  complete  the 
formulation  of  the  problem  we  introduce  the  stream  function  ip,  defined  such  that 


1  dip 
r  dz  ’ 


1  dip 
r  dr  ’ 


(2.25) 


which  makes  the  continuity  equation  automatically  satisfied.  Substitution  of  (2.25)  into  the 
azimuthal  momentum  equation,  (2.17),  results  in 


.  dv  dv  d  ( 1  d  .  ' \ 

- 'f’*“ +  *rTz  =  VTfr 


(2.26) 


Eliminating  it  and  w  in  (2.24)  using  (2.25)  yields 

d  (\ 


Mr,  -  <MM~o  -  (V)  -  *r  [|.  (rA  (4.))  -  |  (r£  (4^ 

=  -r7„'  ^£(rV)dr  ~ 2,/r2  Jo  (;|:(r,’))  ir  +  ’1',T V’" - • 


(2.27) 


7 


Equations  (2.26)  and  (2.27)  are  the  governing  equations  for  the  unknowns  rp  and  v  subject 
to  appropriate  boundary  conditions.  We  note  again  that  (2.27)  is  an  elliptic-integro-partial 
differential  equation. 


3.  Simplified  Model 


We  noted  earlier,  in  a  discussion  of  expression  (2.23),  that  if  the  flow  divergence  a  = 
u/w  >  0,  then  it  follows  that  (dp/dz)T= o  >  ( dp/dz)T-R ,  i.e.,  flow  divergence  in  such  vortex 
flows  can  lead  to  a  less  favorable  pressure  gradient  along  the  axis  compared  to  that  at  the 
edge  of  the  core.  Equation  (2.23)  is  based  on  a  purely  inviscid  argument,  nor  is  anything 
said  about  the  origin  of  the  streamline  divergence.  Within  the  inviscid  assumption  this  may 
most  readily  be  assumed  to  arise  from  an  imposed  adverse  pressure  gradient  at  the  core 
edge,  leading  to  axial  flow  retardation,  and  from  continuity  to  u  >  0.  Could  viscous  effects 
either  cause  or  alter  the  degree  of  axial  retardation,  independent  of  purely  inviscid  effects. 
Grabowski  and  Berger  [14]  conjectured  that  this  could  come  about  as  follows.  Viscous 
dissipation  of  swirl  would  lead  to  a  decrease  in  v  (Eq.  (2.17)),  so  dv/dz  <  0.  It  then  follows 
from  (2.11)  that  there  will  be  greater  adverse  axial  pressure  gradient  and  axial  velocity 
retardation  near  the  axis  than  at  the  core  edge.  This,  according  to  continuity,  (2.15),  will 
lead  to  a  significant  amount  of  outflow  (u  >  0).  This  outflow,  in  turn,  because  of  conservation 
of  angular  momentum,  will  lead  to  a  decrease  in  v  (because  r  increases).  The  combination  of 
outflow,  and  therefore  flow  divergence,  and  diffusion  of  swirl,  leading  to  a  decrease  in  t;  and 
hence  making  dv/dz  <  0  will  combine  to  make  for  a  large  adverse  axial  pressure  gradient 
and  resulting  axial  retardation  near  the  axis.  This  picture  of  the  flow  closely  parallels  that 
of  Hall  [28]. 

The  analysis  presented  in  Section  2  attempts  to  capture  this  physics.  The  resulting  for¬ 
mulation  is,  however,  to>  complicated  for  the  essentially  new  elements,  those  that  distinguish 
this  formulation  from  the  quasi-cylindrical  model,  to  be  evident.  In  an  attempt  to  illustrate 
these,  we  consider  in  this  section  a  simplified  version  of  the  analysis  of  Section  2. 

If  we  evaluate  (2.24)  at  r  =  R,  the  edge  of  the  core,  noting  that  the  viscous  terms  on  the 
left  hand  side  will  be  small  there,  as  well  as  u  dw/dr ,  w  dw/dz,  and  d2w/dz2,  it  reduces  to 


dw  rR  1  /  it  \  d  .  2  2w  [R  1  ( v  \  d  1  d  . 

"'arl"“  =  “/o  ?Ud-r{rv)dr  +  2,/l  ;U  a?  7d-rirv) 


+  |r=0.  (3.1) 


We  anticipate  that  in  situations  where  breakdown  occurs  the  sum  of  the  two  integrals  on  the 
right  hand  side  of  (3.1)  is  negative  since  this  will  lead  to  retardation  on  the  axis  (dw/dz  <  0). 
Because  these  are  integral  terms  they  are  less  sensitive  to  the  details  of  the  velocity  profiles, 
and  as  our  greatest  simplification  we  set  the  sum  of  these  terms  to  a  (negative)  constant 
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(3.2) 


■A,  say.  Equation  (3.1)  then  becomes 


d2w .  dw ,  ,  . 

l/-^r|r=0  -  t«^-|r=0  »  A,  (A>0) 


dz2 


or,  for  simplicity 


dz 

cPw  dw 


-  w—  =  A, 


dz2  dz 

where  w  in  (3.3)  is  understood  to  be  u;|r=o,  the  axial  velocity  on  the  axis. 
First,  we  note  that  for  small  or  negligible  viscosity, 


(3.3) 


with  solution 


dw 

w—~  ~  —A, 
dz 


1  2 

-to  w  —Az  +  const. 

Ld 


(3.4) 


If  w  =  Wc  at  z  =  0,  the  uniform  axial  core  velocity  far  upstream,  then  the  const,  in  the 
above  is  W2 / 2  and 


W  £3  Wr 


For  small  z  this  gives 


to  «  Wc 


(1_  w!z+ 


(3.5) 


(3.6) 


a  linear  decrease  of  axial  velocity  on  the  axis.  This  is  in  qualitative  and  quantitative  agree¬ 
ment  with  Hall’s  calculations  [28]  using  the  quasi-cylindrical  equations. 

3.1.  Exact  solution  of  non-linear  model  equation 

Here  we  present  the  exact  solution  of  non-linear  Eq.  (3.3).  Integration  yields  immediately 

dw  1  ,  /A\  .  .  _ 

Tz  ~  Tvw  =  (t)  z  t  °  (“  =  con8t->-  (3-7) 

This  is  a  Riccati  equation.  It  can  be  transformed  to  a  second-order  1:  ear  ordinary  differential 
equation  by  setting 


w(z)  =  —2u 


y\z) 
y(z )' 


(3.8) 


Equation  (3.7)  becomes 

7?  +  h(7z  +  a)y  =  0'  (3-9) 

which  is  an  Airy  equation.  It  can  be  put  in  standard  Airy  equation  form  by  introducing  the 
new  independent  variable 


v 


WfA 


whereupon  (3.9)  becomes 


(3.11) 


The  general  solution  of  (3.11)  is 


<Py 

SF  -  zv  “  °' 


y(z)  =  cjAi(i)  +  c2Bi(z), 


(3.12) 


where  C\  and  c2  are  arbiiiary  constants  and  Ai  and  Bi  the  Airy  functions.  In  terms  of  the 
original  dependent  variable  w  the  solution  of  (3.3)  is  then 


w(z )  =  (4AI/)1/3 


kAi'jz)  +  Bi'{z) 
kAi(z)  +  Bi(z)  ’ 


(3.13) 


where  the  arbitrary  constants  for  the  general  solution  of  the  second-order  differential  equation 
are  a  and  k. 

As  boundary  or  initial  conditions  for  (3.13)  we  have 


w  =  Wc 


dw 

dz 


,  given 


at  z  =  0. 


(3.14a,  b) 


The  first  of  these  has  already  been  discussed.  As  for  the  second,  far  upstream,  at  z  =  0, 
the  change  in  w  with  2  may  be  small  in  the  absence  of  an  imposed  pressure  gradient,  but 
this  would  not  be  so  if  this  was  not  the  case.  Thus,  for  example,  by  choosing  (dw/dz)0  <  0 
we  can  consider  cases  where  an  adverse  pressure  gradient,  say  as  resulting  from  a  diverging 
channel,  is  imposed  at  the  core  edge. 

Applying  first  (3.14b)  to  (3.7),  we  find  immediately 


a=(%)  -  V’. 

\dz)0  2u 

(3.15) 

Since,  when  z  =  0, 

(  v  y/3 

*z  =  ~\2 T>)  a  =  "°’say’ 

(3.16) 

so,  applying  boundary  condition  (3.14a) 

w  -  a  +  *'(*>) 

c  1  '  kAi(zo)  +  Bi(z0)  ‘ 

(3.17) 

Solving  for  k, 

WcBi(z0 )  -  (4Av)l'3Bi'(z0) 
k  WcAi{z0)  -  ( 4Avy/3Ai'(z0 ) ' 

(3.18) 
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This  completes  the  determination  of  the  two  constants  a  and  k  appearing  in  the  general 
solution,  (3.13),  of  the  non-linear  equation  (3.3).  The  solution  (3.13)  depends  on  the  four 
(dimensional)  parameters  A,u,Wc,  and  ( dw/dz)0 . 

We  note  that  if  [dw/dz)0  =  0,  then 

and  the  solution  simplifies  accordingly. 


3.2.  Asymptotic  solution  for  Reg  — »  oo 


We  consider  here,  using  matched  asymptotic  expansions,  the  asymptotic  solution  of  the 
non-linear  ordinary  differential  equation  (3.3)  for  vanishingly  small  viscosity  i/,  i.e.,  in  the 
limit  Reg  — >  oo. 

We  begin  with  the  outer  expansion,  which  follows  immediately  from  the  first  integral  of 
(3.3),  namely  the  expression  (3.7),  dropping,  to  lowest  order,  the  first  term  on  the  left,  so 

1 


—Az  +  const. 

2 


(3.20) 


If  w  —  Wc  at  z  =  0,  then  the  const.  =  (1/2 )W*,  and  the  leading  order  outer  solution  is 

1/2 


W  ~  Wc 


(  2A  ) 

['~mv  • 


(3.21) 


which  is  exactly  the  earlier  (3.5). 

To  determine  the  scale  of  the  inner  variables,  z  and  w,  we  set 

w  =  Duaw  (a  <  0) 


«  (D<  o) 


where  we  anticipate  both  a,  (3  <  0.  The  differential  equation  (3.3)  becomes 

fpj  _  d-V-2«+0®5  =  AD. 
az ’*  az 

The  two  terms  on  the  left  hand  side  will  be  of  the  same  order  in  v  if 


(3.22a, b ) 


(3.23) 


0  +  a  =  -l. 

If  we  choose  a  and  0  such  that  the  single  term  on  the  right  hand  side  vanishes  in  the  limit 
as  v  — »  0,  then  the  resulting  inner  solution  cannot  be  matched  with  the  outer  solution. 
Therefore  we  must  choose 

—2a  -1-  0  —  0,  or  0  —  2a. 
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It  then  follows  that 


a  =  -1/3,  /?  =  —2/3. 


(3.24) 


Thus  the  scaling  of  the  inner  variables  is 


or  w  =  D  1  vx!zw 


(Wl  \  _2/3  Wc2 

\2A  )  2  A 


(3.25) 


Equation  (3.23)  becomes 


cPw  _  ,  dw  ,  ^ 
-5-  +  D-'w—  =  ^IZ). 
cfz  az 


(3.26) 


This  equation  for  the  inner  region  is  identical,  apart  from  a  change  of  variables,  with  the 
full  equation,  (3.3).  The  following  analysis  therefore  proceeds  along  the  lines  of  that  given 
in  Section  3.1. 

The  first  integral  of  (3.26)  is 


dw  1  ~  x 

—  +  — —  ur  =  Auz  +  a,  (a  =  const.) 
dz  2D 


(3.27) 


This  Ricatti  equation  can  be  converted  to  a  second  order  linear  equation  by  the  substitution 
w  oc  y'/y.  The  solution  of  the  resulting  Airy  equation  is 


where 


W(z)  ~  + 

v  ;  v  ’  kAt(z)  +  Bi(z) 


i=  {2^!1  (AEK  +  a), 


(3.28a) 


(3.286) 


a  and  k  being  arbitrary  constants,  to  be  determined  by  matching  of  the  inner  and  outer 
solutions. 

To  match  the  outer  solution  (3.21)  and  the  inner  solution  (3.28),  we  begin  by  rewriting 
the  outer  solution  (3.21)  in  inner  variables  and  expanding  for  v  — *  0,  with  z  held  fixed, 
obtaining 

~  2Azp2/3.  (3.29) 


Rewritten  in  outer  variables  the  inner  variable  z  is 


4-2/3 

2V3  D 


(S-h 


-fa. 


(3.30) 
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For  reasons  discussed  earlier  A  >  0,  and  we  can  assume  that  D,  a  scale  factor  for  w,  is  also 
positive.  It  then  follows  that  in  the  limit  v  — ►  0,  with  z  held  fixed,  that 


W2 

2  — »  oo  if  —  z  >  0, 
2A 

W2 

z  — *  — oo  if  —  2  <  0. 
2A 


(3.31) 


If  the  second  of  these  alternatives  is  the  case,  Ai(z)  and  Bi(z)  oscillate  and  the  inner  solution 
(3.28)  cannot  be  matched  to  the  outer  solution,  (3.29).  Therefore  we  assume  the  first  case, 
that  z  — ►  oo  as  v  — >  0.  For  asymptotically  large  positive  arguments  Ai  and  Bi  have  the 
following  asymptotic  expansions 


'  Ai(S)  „  _i_r,/4e-2i3/2/3,  f  Bi(z)  ~  -Lz-'/'e2*312'3, 

2  y/n  y/TT 

as  z  — *  oo  <  < 

Ai'(z)  ~  -ZliVV253'2/3,  Bi\z)  ~  4=^1/4e2i3/2/3- 

2-y/ 7 T  ^  y/  7T 

It  then  follows  that  in  this  (outer)  limit  that  (3.28a)  behaves  asymptotically  as 

w(z)  ~  (4 A)1/3Dz1'2, 
so,  according  to  (3.25)  and  (3.28b), 

Winner  ~  (4A)1/3U1/3  Zl/2 


(3.32) 


(3.33) 


(4  A)1^3!/1 


(2A)-1'3 


{ADz  +  a) 


WLer  ~  ( 2Az  +  )  V213. 


(3.34) 


The  variable  part  of  the  inner  and  outer  solutions,  (3.29)  and  (3.34)  agree  exactly. 


3.3.  Exact  solution  in  the  limit  of  vanishingly  small  viscosity 

We  consider  the  limit  of  the  exact  solution  given  in  Section  3.1  as  v  — ►  0.  In  this  limit, 
from  (3.15)  and  (3.16), 


>  /  1/  \  J/3  /  1 


*-(*)• -Gad  — 


(3.35a) 


(3.356) 
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From  (3.18),  using  the  asymptotic  expressions  (3.32), 


lim  k  =  0. 

v->0 


(3.35c) 


From  (3.10) 


lim^o  z 


_(J_\'I3(A  l 


\2AV 


\l>  2v  J 


\(2A2u2)~1I3(W2  -  2 Az) 


[  +oo,  if  W2  >  2Az 


(3.36a) 


(3.366) 


[  -oo,  if  Wc2  <  2 Az. 

We  assume  that  W2  >  2Az ,  so  z  — >■  -foo  as  v  — >  0.  It  then  follows  from  (3.13),  using  (3.35c), 
(3.32)  and  (3.36a),  that 

1 1/2 


ri  / 

limtc(i)  — >  (iAuY^z1^2  =  (\Av)1^  -(2A2u2)~1^3  (w2  —  2Az)  , 


(3.37) 


(W2  -  2Azy'2, 


which  is  identical  with  the  outer  solution,  (3.21),  presented  in  Section  3.2  as  part  of  the 
asymptotic  solution  for  Reg  —*  oo,  and  with  the  limiting  small  viscosity  solution  given 
earlier  (Eq.  3.5). 


3.4.  Non-dimensionalization 


We  begin  our  non-dimensionalization  of  (3.3)  by  first  considering  the  quantity  on  the 
right  hand  side,  A ,  which  stands  for  the  (negative  of  the)  two  integral  terms  on  the  right 
hand  side  of  (3.1).  We  consider  each  of  these  in  turn.  The  order  of  magnitude  of  the  first  of 
these  can  be  estimated  as  follows: 


— 


(3.38) 


where  rv  =  k  ~  T,  the  characteristic  circulation  of  the  vortex,  8  the  characteristic  core 
radius  (~  R)1  and  u/w  ~  a,  a  characteristic  flow  divergence. 

The  magnitude  of  the  second  integral  on  the  right  hand  side  of  (3.1)  is  similarly  estimated: 


2"  /:  j  (^  i  ;  ih  dr  ~ 2"  ©  © 6 = 2"  © r-  <3-39> 


where  v/w  ~  s,  a  characteristic  swirl  velocity  ratio.  (Both  integrals  have  dimensions  of 
LT~ 2.)  Unlike  (3.38)  which  is  generally  negative  as  the  vortex  approaches  breakdown,  the 
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sign  of  the  second  integral  is  not  obvious.  To  allow  us  to  consider  the  effect  of  the  integral 
being  negative  (driving  the  vortex  to  breakdown)  or  positive,  we  shall  attach  to  it  an  arbitrary 
sign,  denoted  by  sgn.  Then 


^r2  +  2sgm,(£)  r> 


and  the  differential  equation  (3.3)  becomes 

<Pw  dw  ( ct 


diw  dw  /  a\  ~  „  /  s  \ 

"-d?  -  = [p) r  -2^yr- 


(3.40) 


(3.41) 


We  now  non-dimensionalize  (3.41)  by  defining  non-dimensional  variables  and  parameters 


z  _  w 
2  =  w  =  777 


r  D  wcs 

7_  SWZ  Re‘~  v  ' 


(3.42) 


with  which  definitions  (3.41)  becomes 


1  <Pw  _ dw  2  2 

rTsI?  -  W~E  =  07  +  Ws  sgni7' 


(3.43) 


The  solution  of  (3.43)  depends  on  the  four  parameters  Re$ ,  the  core  Reynolds  number, 
a(~  u/w ),  the  flow  divergence,  s(~  v/w),  the  swirl  velocity  ratio,  and  7(=  r/£Wc),  the 
non-dimensional  flow  circulation. 


4.  Calculations  and  Results 

We  have  calculated  a  number  of  solutions  of  Eq.  (3.43),  by  evaluation  of  the  exact 
solution  given  in  Section  3.1,  in  order  to  exhibit  the  dependence  of  the  solutions  on  the 
four  parameters  Res ,  a,  s,  and  7  which  appear  in  the  equation,  and  the  value  of  the  initial 
velocity  gradient  (dw/dz) 0.  Figures  1  and  2  show  the  non-dimensional  axial  velocity  variation 
tu(z)  with  Re$  (for  convenience  the  bars  over  the  non-dimensional  variables  introduced  in 
Section  3.4  have  been  removed)  for  the  same  values  of  (dw/dz) 0(=  —0.1)  and  o(=  0.005) 
but  different  values  of  s(=  0.5  and  1.0)  and  7(=  1.0  and  1.5).  (We  anticipate  that  negative 
values  of  (dw/dz) 0,  corresponding  to  an  adverse  axial  pressure  gradient,  are  most  likely  to 
lead  to  breakdown;  this  is  borne  out  by  sample  calculations  of  such  cases.  Also,  in  all  the 
calculations  reported  here  the  sgn  in  (3.40)  is  taken  to  be  positive,  implying  that  the  integral 
in  (3.39)  is  negative,  as  is  that  in  (3.38),  increasing  the  tendency  of  the  vortex  to  undergo 
breakdown.)  We  note  that  for  both  cases  in  Figs.  1  and  2  the  vortex  undergoes  breakdown, 
i.e.,  the  axial  velocity  becomes  zero,  for  all  Res  in  the  range  50  <  Res  <  1000,  except 
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possibly  for  Reg  =  1000.  Moreover,  the  point  on  the  axis  where  stagnation  or  breakdown 
occurs  moves  rapidly  upstream  as  Reg  increases,  reaching  asymptotically  a  position  close  to 
the  initial  station.  Finally  we  note  that  the  curves  in  Figs.  1  and  2  are  remarkably  similar  for 
the  same  Reg ,  suggesting  that  the  results  are  insensitive  to  the  values  of  s  and  7.  Calculations 
confirm  this  for  other  values  of  the  parameters  as  well.  Figure  3  shows  the  axial  velocity 
for  the  same  Reg  range  and  the  same  s  and  7  as  in  Fig.  2  but  a  smaller  value  of  a  and  a 
larger  negative  initial  axial  velocity  gradient.  We  note  in  comparing  Fig.  3  to  2  the  much 
more  rapid  decrease  in  axial  velocity  as  Reg  increases.  This  is  primarily  a  consequence  of 
the  much  larger  negative  value  of  ( dw/dz)o  in  Fig.  3.  (That  the  change  in  a  has  only  a  small 
effect  can  be  seen,  for  example,  from  Fig.  4  for  the  same  values  of  all  the  parameters  as  in 
Fig.  3  except  for  a  which  is  ten  times  as  large,  0.01  as  compared  to  0.001.  There  is  very 
little  discernible  difference  between  the  curves  in  Figs.  3  and  4.)  Figure  5  shows  the  same 
curves  for  the  same  values  of  the  parameters  as  in  Fig.  2  except  that  for  Fig.  5  |(dio/dz)0|  is 
a  tenth  as  large.  The  consequence  of  this  is  that  compared  to  Fig.  2  the  curves  are  displaced 
further  downstream,  and  even  more  significantly,  for  the  smaller  Reg  the  axial  flow  does  not 
stagnate  but  remains  relatively  constant  and  then  increases  with  axial  distance!  Again  this 
illustrates  the  importance  of  the  value  of  the  initial  axial  velocity  gradient. 

Because  of  the  particular  significance  of  the  value  of  the  axial  location  where  the  axial 
velocity  becomes  zero,  denoting  the  location  of  breakdown,  we  have  used  the  method  of  con¬ 
tinuation  [29]  to  determine  the  variation  of  this  point  as  a  function  of  the  various  parameters. 
The  results  above,  as  well  as  others  unreported  here,  demonstrate  that  the  point  at  which 
the  axial  velocity  becomes  zero,  denoted  by  zv_ varies  primarily  with  the  parameters  Reg 
and  ( dw/dz)0 ■  Figure  6  shows  the  variation  of  zv,b.  with  Reg  for  four  values  of  —(dw/dz)0\ 
the  other  parameters  have  the  values  a  =  0.01,  7  =  1,  and  s  =  0.5.  We  see  that  for  all  these 
values  of  the  initial  gradient  the  location  of  axial  stagnation  moves  rapidly  upstream,  and 
that  while  the  exact  values  of  zv, asymptote  for  larger  Reg  there  are  significant  differences 
for  smaller  Reg.  Figure  7  shows,  for  Re  =  500,  a  =  0.015,  7  =  1,  and  s  =  0.5,  the  variation 
of  z„.6.  with  —(dw/dz) 0.  We  note  the  very  rapid  upstream  migration  of  the  location  of  axial 
breakdown  as  the  initial  axial  velocity  gradient  becomes  more  negative,  i.e.,  larger  initial 
adverse  pressure  gradient. 

5.  Conclusions 

The  principal  aims  of  this  paper  were  to  show  the  origin  of  the  pronounced  tendency  of  the 
location  of  vortex  breakdown  to  migrate  upstream,  an  effect  characteristic  of  most  numerical 
simulations  and  noted  in  some  experiments.  Based  on  an  order-of-magnitude  analysis  of  the 
viscous  Navier-Stokes  equations  a  modified  set  of  governing  differential  equations  for  swirling 
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flows  has  been  derived.  A  theoretical  analysis  of  these  equations  shows  that  unlike  the 
parabolic  quasi-cylindrical  equations,  which  they  resemble,  these  equations  have  an  elliptic 
character.  This  analysis  also  exhibits  how  the  location  of  axial  stagnation,  interpreted  here 
as  a  signal  of  breakdown,  depends  explicitly  on  a  number  of  parameters:  the  core  Reynolds 
number,  the  flow  divergence,  the  swirl  velocity  ratio,  and  the  strength  of  the  vortex.  A 
simplified  model  for  the  axial  velocity  variation  is  presented  and  solved  exactly.  Numerical 
results  obtained  from  this  solution  for  a  range  of  the  governing  parameters  show  that  there 
is  a  very  pronounced  tendency  of  the  breakdown  location  to  move  upstream  as  the  core 
Reynolds  number  increases,  or  the  initial  adverse  pressure  gradient  increases,  in  qualitative 
accord  with  most  numerical  and  experimental  simulations. 
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Figure  1.  Variation  of  axial  velocity  with  axial  distance  for  different  core  Reynolds  numbers; 
parameters:  a  =  0.005,  7  =  1.5,  s  =  1.0;  initial  axial  velocity  gradient,  ( dw/dz)0  =  —0.1. 
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Figure  2.  Variation  of  axial  velocity  with  axial  distance  for  different  core  Reynolds  numbers; 
parameters:  a  =  0.005,  7=  1,  »  =  0.5;  initial  axial  velocity  gradient,  ( dw/dz)o  =  -0.1. 
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Figure  3.  Variation  of  axial  velocity  with  axial  distance  for  different  core  Reynolds  numbers; 
parameters:  a  =  0.001,  7  =  1,  a  =  0.5;  initial  axial  velocity  gradient,  ( dw/dz)0  =  -0.15. 


2 

Figure  4.  Variation  of  axial  velocity  with  axial  distance  for  different  core  Reynolds  numbers; 
parameters:  o  =  0.01,  7=1,  a  =  0.5;  initial  axial  velocity  gradient,  (dw/dz)o  =  -0.15. 
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Figure  5.  Variation  of  axial  velocity  with  axial  distance  for  different  core  Reynolds  numbers; 
parameters:  a  =  0.005,  7  =  1,  s  =  0.5;  initial  axial  velocity  gradient,  ( dw/dz)0  =  —0.01. 
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Figure  6.  Location  of  the  axial  vortex  breakdown  point  as  a  function  of  core  Reynolds 
number  for  different  initial  axial  velocity  gradients;  parameters:  a  =  0.01,  7  =  1,  3  =  0.5. 
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Figure  7.  Location  of  the  axial  vortex  breakdown  point  as  a  function  of  the  initial  (negative) 
axial  velocity  gradient;  parameters:  a  =  0.015,  7  =  1,  s  =  0.5,  Res  =  500. 
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